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ABSTRACT 

A multi -input/output random vibration control algorithm was developed based on 
system identification concepts derived from random vibration spectral analysis 
theory. The unique features of the algorithm are: 1) the number of input exci tors 

and the number of output control responses need not be identical, 2) the system 
inverse response matrix is obtained directly from the input/output spectral matrix, 
and 3) the system inverse response matrix is updated every control loop cycle to 
accommodate system amplitude nonlinearities. A laboratory demonstration case of two 
inputs with three outputs is presented to demonstrate the system capabilities. 

INTRODUCTION 

All products are subjected to vibration during their lives, either as a result 
of the manufacturing process, transportation or in-service environments. The 
capability of the product to withstand this vibration environment and to provide the 
desired reliability becomes a major design factor in high-technology applications 
where designs are constrained by severe environments, limited space envelopes, 
weiqht constraints and cost. These designs are often beyond the state-of-the-art in 
analysis procedures, and one is forced to simulate the product s environment to 
determine its reliability. To accurately estimate in-service reliability, it then 
becomes necessary to develop a simulation program which accurately reproduces fie 
environments. 

In the most general sense, a field vibration environment will consist of motion 
with six degrees of freedom. Presently available experimental hardware and control 
systems are limited in their capability to simulate this environment. Most tests 
are performed using uniaxial random excitation and repeating the excitation along 
three mutually perpendicular axes. The primary concern is that the mean-time- 
between-fai lure (MTBF) rates determined from uniaxial tests do not represent in- 
service MTBF rates. One problem is determining the duration of testing along each 
axis of a uniaxial test and the respective test level adjustments to account for 
physical coupling actually experienced in the field. Because of these 

uncertainties, the MTBF obtained from uniaxial test can only be considered an 
estimate of the in-service MTBF. For aerospace product development where product 
weight can have considerable influence on overall system performance, the need for 
reliable design evaluation is critical. A secondary consideration is the cost of 
three test setups to complete the simulation cycle in a uniaxial test. 

To satisfy these shortcomings, it will be necessary to develop test facilities 
that simulate the field environments using multi-axis excitation. In the words of 
MIL-STD-810 (Ref. 1), the latest revision of the military's environmental test 
procedures, the test should be 'tailored' to measured field vibration. This will 
require multi-axis vibration systems with control of each axis and their 
interaction. 
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The objectives of the present study were to develop a multiple-input/output 
random vibration control algorithm based on system identification concepts derived 
from random vibration spectral analysis theory and to evaluate the usefulness of the 
algorithm via a laboratory-based demonstration. 

A number of researchers have described in various detail the historical 
development of multi-axis/actuator control systems [Refs. 2-4]. Fisher at Lawrence 
Livermore Laboratory was one of the first to integrate a control theory and the 
required hardware to demonstrate system capabilities and limitations to perform 
multi-shaker vibration simulation [5]. The digital control system utilized a cross- 
coupling compensation matrix for signal control in the dc to 500 Hz frequency 
range. The major limitations of the system was the orthogonal axis excitation and 
the restriction that the number of response points must equal the number of 
excitation points. 

Several multi-axis shaker systems are in present use; however, a number of them 
employ independent actuator control without cross-axis compensation. One exception 
is the. facilities at Sandia National Laboratory [6-8] which has two-shaker random 
vibration control system with software capable of controlling up to four shakers. 
The control algorithm is based on determination of the inverse system frequency 
response matrix, which requires that the number of inputs and outputs must be 
identical. No capability exists for updating the system transfer function matrix 
directly during the control process, however, an incremental correction matrix is 
used to improve control. 

There are several desired characteristics of a random vibration control 
algorithm that are believed to be important in simulating field vibration 
environments, they are: 1) Provide true stationary random excitation with near 
Gaussian amplitude distribution and signal crest factor (peak to r.m.s.) control; 
2) Full control of cross-coupling between all output responses; 3) Number of input 
signals and number of output responses need not necessarily be the same; 4) 
Maintain loop equalization time to within a small fraction of total test time (say, 
less than 10%); and 5) Provide sufficient dynamic range and frequency bandwidth to 
meet simulation tolerances. Only when multi-axis testing becomes a reality can the 
degree of importance of each of these characteristics be determined. 

In the sections to follow, a random vibration control algorithm is described 
that removes many of the limitations of existing algorithms, the software/hardware 
implementation of the algorithm is described, and an example of test results using 
the system are given to point out its capabilities and limitations. Several 
additional examples are given in Ref. 9. 

CONTROL ALGORITHM FORMULATION 

The development of the random vibration control algorithm views the simulation 
problem as one of system identification. The approach uses the spectral analysis 
techniques discussed in Ref. 10. For the purposes of discussion, consider the two- 
input/three-output simulation problem shown in Fig. 1. Physically, the "system" 
represents the dynamic response between voltage signals input into the vibration 
exciters, denoted as x-^(t) and x 2 (t), and the response accelerations of an item 
attached to the exciters, denoted as y^t), y 2 (t) and y 3 (t). The finite period (T) 

Fourier transform of the kth sample record is denoted as X^(f,T), i = 1, 2; 

and Yj(f,T), j = 1, 2, 3, respectively. The signal Fourier transforms are used to 

construct one-sided input/output cross-spectral densities as 
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( 1 ) 


G xj( f ) ■ jl!f e 

with input auto-spectra, 

G xx (f) = t!»T E * x 5( f ’ T >! (2) 

The output auto-spectra G yy(f) is constructed in a similar way. The operator E is 

the mathematical expectation, and when averaged over N sample records of length T, 
the one-sided spectral density functions take the general form 

G xj (f) = FT E [X i (f ’ T) * * Y j (f * T)1 (3) 

where the * denotes the complex conjugate operator. The signals are assumed to be 
stationary random and, therefore, a finite number of samples, N, may be used to 
estimate the desired cross- and auto-spectra. 

The system frequency response function, [Hy X (f)J, is defined by 

{Y(f)> = [Hy X (f ) 3 {X ( f ) } (4) 

where the brackets { } denote the vector of response or input Fourier spectra. The 
inverse frequency response function, [I X y(f)], is defined as 

(X (f ) } = l I xy ( f ) ] ( Y ( f ) } (5) 

where the two system matrices are generally considered to be the inverse of each 
other. Such a concept is only valid when the number of input channels equals the 
number of output channels, which is not the general case considered in present 
development. To estimate [ I xy ( f ) ] » we use the concept of mathematical expectation 

applied to Eq. 5 after post multiplying by the conjugate transpose of { Y(f ) } , which, 
after some algebra, results in 

[ V f)) = [G xy( f > 1[G yy (f)rl - (6) 

Thus, the system inverse response matrix is the product of the system input/output 
cross-spectral matrix and the inverse of the system output spectral matrix. It is 
important to note that when the system outputs are fully correlated, the output 
spectral matrix is rank deficient and the usual inverse is not defined. When this 
occurs one could drop control of one of the responses, combine two responses in an 
average or least squares way or define a pseudo inverse. The present algorithm 
employs a pseudo inverse as described by Greville [11]. The pseudo inverse routine 
produces the true inverse for rank sufficient matrices while producing unique 
inverses for rank deficient matrices. Due to extraneous noise at the output of most 
physical systems the pseudo inverse is often not necessary. 

The relationship between the input auto-spectra and output auto-spectra is 
derived in a similar way as 

l G x X (0] = [I* y (OHG yy (f)][I xy (f)i T . (7) 
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Given a required output auto-spectra [Ry y (f)3> which is generally an envelope 

of multiple field measurements, the required input auto-spectra is defined by 
Eq. (7) as 

i G xx (f)i ■ [i xy (f)HR yy (f)][i x:y (n] T . (8) 

The control problem becomes that of defining [X (f ) } such that the [G xx (f ) ] of 

Eq. (8) is realized. If the input signals were statistically independent random 
signals, there would be no output cross-coupling control other than that fixed by 
the physical exciters and test item arrangement. Thus, it becomes necessary to 
introduce cross-coupling control. This is accomplished by considering a coupling 
matrix [A(f ) ] such that 


(X(f)} = [A(f ) ] {$(f )} , (9) 

where the [ A(f ) ] matrix is lower triangular (upper triangular elements are all zero) 
and the input {S (f ) } vectors are mutually independent white noise sources (unit 
r.m.s. amplitude with random phase). Upon post multiplication of Eq. (9) by the 
conjugate transpose of (X(f)], and taking the mathematical expectation, there 

results an expression for the coupling elements in [A(f)], 

[A(f)][G ss (f)][A*(f)l T = [G xx (f)]. (10) 

Since the input {S(f ) } vectors are mutually independent white noise, their spectral 
matrix [G s _(f)j is a unit diagonal matrix. A solution for the coupling elements in 
[A] is given on page 262 of Ref. 10. 

A flow chart is given in Fig. 2 which summarizes the major tasks of the control 
algorithm. To insure test item safety during control signal shaping, an incremental 
buildup to the full level drive has been introduced as denoted by the constant a and 

increment A a. As can be seen in Fig. 2, system amplitude nonlinearities are 

compensated via iteration by updating the system inverse response matrix [I xy (f)l 

every control loop. Maximum control is often achieved within two to three loop 
cycles as will be discussed in the results section. 

Overall spectral error E yy is generated for each output by computing the mean 

and standard deviation of the difference in logarithmic levels between the desired 
level R yy and the measured level G yy across the frequency spectrum. The maximum 

high and low values are also monitored. If control appears to be adequate, with 
less than a 2 to 3 dB mean variation, or improvement does not seem to be occurring, 
then individual spectrum plots are made to visually inspect the overall spectral 
matching achieved. The mean and standard deviation error estimates are useful for 
reducing the time consuming graphics data display during loop build up to the final 
drive level. Final judgment on acceptable drive levels are always made based on 

graphical inspection of the data. 

SOFTWARE/HARDWARE IMPLEMENTATION 

The hardware configuration used to implement the control algorithm was 
configured from two main stand alone processors. The drive signal generation and 
overall loop control is carried out via a 32-bit virtual memory laboratory 
computer. The laboratory computer generates the drive signals via a digital-to- 
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analog (DTA) interface and supplies control signals to a multi-channel fast Fourier 
transform (FFT) analyzer. The multi-channel FFT analyzer is capable of digitizing, 
Fourier transforming, and sample averaging up to 8 channels of data while 
maintaining time correlation. The input drive signals are sampled simultaneously 
with the output response signals so as to maintain correlation of all input and 
output signals. The various sample-averaged cross-spectral densities required by 
the control algorithm are transferred from the FFT analyzer to the laboratory 
computer via direct memory access (DMA). The transfer of one channel pair of auto- 
and cross-complex spectra (400 line) requires approximately 50 milliseconds. A two- 
input/four-output system requires transfer of 14 channel pairs of information, thus 
requiring approximately 0.7 seconds. The spectra transfer generally requires far 
less time than the sample averaging time required for a confident data sample when 
considering random error. 

Drive Signal Generation 

As previously stated, the approach taken to generate the required input .drive 
signals is to consider the shaker inputs to be a combination of statistically 
independent random, stationary Gaussian signals {S^}, i=l,2,3 ...M. In which case, 

the required drive signals {Xj}, 3=1, 2, 3 ...N can be related to the ideal signals 

{S} via a lower triangular coupling matrix [A]. The initial, drive signals are 
generated first in terms of their spectral content; namely, unit amplitudes exist 
from the lowest frequency, f^_, to the highest frequency, f^, of interest for output 

control. Using a uniformly distributed random number generator, the phase of each 
spectral component is generated with random variation from component to component. 
To be consistent with the FFT analyzer, a maximum of 400 spectral lines exist 
between zero frequency and the analyzer cutoff frequency. The analyzer cutoff 
frequencies, f^, are predefined and selected such that f|^ < f^. As such, the 

Fourier spectrum of the drive signals initially consist of unit amplitude random 
phase components in the analysis range f L to f H . Finite length time histories of 

the initial signals are then obtained via the inverse Fourier transform of each 
independent signal. 

In order to capture a number of independent random sample averages of the 
system response, sequences of the time histories must be continuously generated. If 
the same time series for a particular channel is repeatedly generated, the output 
wave form becomes pseudo-random and not acceptable for environmental simulation. To 
form acceptable drive signals, sequences of the same signal were randomly delayed by 
approximately 25 percent of the signal period and summed for output. It was. found 
that if the signal time delay were fixed and not randomized, the resulting signals 
showed marked periodic components with unacceptable coherence between "independent" 
drive signals, as is shown in Fig. 3. However, when introducing the randomized time 
delays, acceptable signals such as shown in Fig. 4 resulted. It was found that the 
random phase relationship between the original independent drive signals {S} could 
only be maintained when introducing the coupling matrix in the following partitioned 
form (as for a two-input system); 



143 



Time histories are generated for the X 11# X 21 and X 22 signals individually where X^ 
and X 21 maintain identical random phase shifts and random time lags as specified by 
the original generation of S^, and X 22 that of $ 2 . The desired time history drive 
signals corresponding to the specified coupling (the A matrix) are then 

*l(t) = x u (t) 

and (12) 

*2( t ) = x 21 (t) + x 22 (t). 

Additional important characteristics of the drive signals are that they be 
stationary and that the signal crest factor can be controlled. Stationarity is 
measured by the variation in signal mean and r.m.s. level as time progresses. For 
the present algorithm, the zero mean and r.m.s. variation are less than 1 percent. 
Control over peak-to-r.m.s. level is accomplished by introducing signal peak 
folding. All data points exceeding the desired peak value are folded back into the 
signal, not clipped at the peak value. Thus, signal peaks when folded become local 
valleys which tend to raise the signal r.m.s. level. Crest factor control in the 
range from 2.0 to 4.5 is readily accomplished. 

Output Spectra Specification 


In order to accommodate various levels of output control specification which 
presently exist for equipment under qualification and provide a means to evaluate 
the importance of all elements in a specification, the following specifications for 
the Ryy spectral matrix were developed. 


1) 


2 ) 


Full Ryy Specification: where auto-spectra are specified by spectrum 

frequency and magnitude break point pairs and cross-spectra are specified 
by spectrum frequency, magnitude, and phase break point values. 

Auto-Spectra and Coherence Specification: where auto-spectra are specified 

by spectrum frequency and magnitude break point pairs and cross-spectra are 
specified by frequency and coherence break point pairs. Herein, the cross- 
spectra magnitudes are generated from the definition of spectral coherence. 


Y 


2 

i j 


G. . 
II 


G..G.. 
n JJ 


(13) 


and the phase of R yy is taken as that of the measured G yy . 

3) Auto-Spectra Only: where auto-spectra are specified by frequency and 

magnitude break point pairs. The cross-spectra magnitude and phase are 
generated from the measured G yy in proportion to the magnitude of the 

specified R yy via: 



( R ii • V 


(G.. • G..) 



(14) 


For existing test specifications such as MIL-STD-810C/D, only the output auto- 
spectra may be known. If field measurements were made without specification of 
cross-spectra, one may be much better off to accept the cross-axis coupling of the 
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test setup rather than attempting to invoke a physically unrealizable specification 
for which control would be impossible or which may drastically alter the service 
life of the test item. 

TYPICAL RESULTS 

Several types of exciters, exciter arrangements, and test item configurations 
were used to evaluate the proposed control algorithm and its limitations as reported 
in Ref. 9. The algorithm evaluation process will, of course, continue as the system 
is employed for equipment qualification purposes, and the author is confident that 
much is to be learned as opportunities arise for its use. The following example 
illustrates the general capabilities of the control algorithm. 

Multi -Resonant Beam 

In order to study the control of sharp structural resonances a small cantilever 
beam was driven using two small voice coil shakers. Three accelerometers were used 
to record the beam response. The two-input/three-output system is schematically 
shown in Fig. 5. 

The cantilever beam exhibited damping of less than 0.3% critical, which is 
uncharacteristic of built-up structures in which damping levels of 2 to 5% are more 
representative. Thus, the low damped beam provides a harsh test of the control 
algorithm to control resonances. The initial drive signals to the exciters are 
shown in Fig. 6 as being independent random signals in the frequency range from 5 to 
250 Hz. The initial beam response to the drive signals and the desired output 
magnitude specifications are shown in Fig. 7. Due to the sharp resonances in the 
beam, exhibiting distinct phase between outputs, only the cross-axis coherence was 
specified since the excitors were of marginal power and, thus, not likely powerful 
enough to change system phasing. The coherence between output channels 1 and 2 was 
set at 0.50, as was that between channels 2 and 3. The coherence between channels 1 
and 3 was set to 0.20. Control after three loop iterations resulted in the output 
spectra shown in Fig. 8. As can be seen, the primary beam resonance was well- 
controlled, as was the second beam resonance at outputs 2 and 3. Since output 1 was 
at an anti-node of the second beam resonance, sensitivity to the second mode was 
high at that location as is shown in Fig. 8a. The second beam resonance occurs 
around 90 Hz and with a critical damping ratio of 0.003, its half-power point 
bandwidth would be approximately 0.54 Hz. Thus, with a control bandwith of 0.625 Hz 
(250 Hz/400 lines), reasonably good control of the mode was realized. Control of 
the cross-spectra was reasonably good, however, not as good as the auto-spectra 
shown in Fig. 8. 

Aircraft External Store Vibration Simulation 

The aircraft external store used to evaluate a two-input/four-output 
control configuration is shown in Fig. 9. The external store is a dummy test item 
with an overall length of 193.3 cm (76.12 inches). Major bulkheads are at Missile 
Stations (MS) 30.38, 55.06 and 66.125. The bulkheads provide for exciter attachment 
and free-free missile suspension via elastic bungy cords (2 Hz vertical). The test 
item is fitted with add-on damping material to raise modal damping values to 0.02 
minimum. A dynamic modal analysis of the test item in the frequency range from 10 
to 500 Hz revealed an abundance of normal mode resonances as shown in Table 1. It 
was felt that this multi-modal structure would provide a good evaluation of the 
proposed control algorithm. 
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Input excitation to the test item was provided by two electrodynamic modal 
shakers attached at MS 30.38 and MS 66.125. The shakers were aligned normal to the 
missile longitudinal axis and 45 degrees to the vertical (see Fig. 9). The 
acceleration response control points are located at MS 55.06 and MS 30.38 with one 
vertical and one horizontal transducer at each location. Response channels 1 and 4 
are the horizontal responses and channels 2 and 3 the vertical responses. 
Field data were not available for the test item and, therefore, the control 
specification was a flat response in the frequency range from 20 to 500 Hz on each 

4 2 

output with magnitude 10 g /Hz. The bandwidth of resolution was set at 1.25 Hz 
which should provide marginal resolution to control a 2 % damped 100 Hz mode (half 
power point bandwidth of 4.0 Hz). 

The output response specification for channels 1 and 3, along with the recorded 
responses after the initial broadband excitation and after three control loop 
cycles, are shown in Fig. 10. As can be seen by the data given in Fig. 10a, 

channel 1 response was initially dominated by the fuselage first and second 
horizontal modes; nevertheless, the control algorithm performed reasonably well. 
Similar trends were seen in the other responses as well as is shown in Fig. 10b for 
the channel 3 response. Close examination of the hay-stack response in the area of 
460 Hz in channel 3 response shows the incompatibility of the control specification 
relative to the location of the excitors and response accelerometers. In this case, 
it was impossible to drive response 4 (not shown) in the horizontal direction to the 
desired level while maintaining vertical response to within limits in the frequency 
region of the third fuselage vertical bending mode at 462.5 Hz. For the most part, 
the control algorithm performed well, except at these incompatible areas. It is 
expected that measured field data would not reflect these incompatibilities if the 
exciters were placed at the store-to-aircraft attachment points and field data were 
recorded at the specified output locations. 

Throughout the evaluations discussed above, 20 sample averages were used to 
obtain the results shown. In many cases, 5-10 averages would have been adequate 
during signal development. A higher number of averages, on the order of 200, are 
appropriate for random data analyses of the final drive signals to establish 
statistical confidence [10]. Control loop time is highly dependent on the frequency 
range of interest and the graphical display of results. Twenty sample averages at a 
maximum frequency of 62.5 Hz requires 128 seconds; while for an analysis out to 
2000 Hz, the sample average time is only 4 seconds. Thus, in the former case, 10 
seconds for control analysis time is insignificant while in the latter case it is 
measurable. In either case, the time to record graphical data to a hardcopy device 
is significant. Such intermediate output records can, of course, be avoided after 
developing confidence in the algorithm. 

CONCLUSIONS AND RECOMMENDATIONS 

Based on the limited experience gained with the control algorithm thus far, the 
following conclusions and observations are made: 

1) It appears that the control algorithm, as implemented, will allow different 
numbers of input and output channels. This feature can be useful in 
various practical test situations, including assurance of uniform 
excitation for a test item. 

2) The algorithm updates the system inverse response matrix every control loop 
cycle which accommodates system amplitude nonlinearities and allows 
compensation for system degradation during life time simulation. 
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3) Full control of cross-coupling between all output responses can be invoked 
with the present algorithm. However, it appears that the extent to which 
control is achieved is highly dependent on the physical arrangement of the 
exciters. 

4) Analysis bandwidth is an important parameter in controlling test item 
resonances. The 400-line analysis used in the present system appears to be 
marginal, with 800 to 1000-line analysis being preferred. 

5) The dynamic range of the present system is set by the FFT analyzer, that 
being the 12-bit analog-to-digital converters. It is recommended that 16- 
bit resolution be used to insure sufficient dynamic range especially for 
random signal control where crest factors of 3 to 4.5 are required. 

SYMBOLS 


Values are given in both SI and U.S. Customary Units. The measurements and 
calculations were made in U.S. Customary Units. 
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A(t) 
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X) 

e$f> 

G y| (f) 

H yx (f) 

I xy (f) 

N 

S(t) 

T 
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Table 1. External Store Natural Frequencies 


Mode No. 

Frequency 

(Hz) 


Mode Description 

1 

15.0 

1st 

Symmetric Wing Vertical Bending 

2 

52.5 


Symmetric Wing Torsion 

3 

57.5 


Coupled Tail /Wing Horizontal 

4 

97.5 

1st 

Fuselage Vertical Bending 

5 

102.5 

1st 

Fuselage Vertical Bending 

6 

107.5 

1st 

Fuselage Horizontal Bending 

7 

120.0 

2nd 

Symmetric Wing Vertical Bending 

8 

202.5 

2nd 

Fuselage Horizontal Bending 

9 

227.5 

2nd 

Fuselage Vertical Beding 

10 

390.0 


Wing Panel Mode 

11 

410.0 


Wing Panel Mode 

12 

435.0 

2nd 

Fuselage Horizontal Bending 

13 

462.5 

3rd 

Fuselage Vertical Bending + 2nd Wing 
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Coherence Power Spectral Density 

Volts 2 /Hz 



Frequency -Hz 
a) Auto Spectra Channel 2 



Frequency -Hz 

b) Coherence Between Channels 1 and 2 


Figure 3. Typical Spectral Averaged Drive Signals, Fixed Time Delay 
During Time Sequence Generation 
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Frequency -Hz 
a) Auto Spectra Channel 2 
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Figure 4. Typical Spectral Averaged Drive Signals, Random Time Delay 
During Time Sequence Generation 




Coherence Power Spectral Densit 


Accelerometer Outputs 



Force Coil Inputs 

Figure 5. Schematic of the Two-Input/Three-Output System 



a) PSD -Channel 1 or Channel 2 



Frequency -Hz 

b) Coherence Channel 2 to Channel 1 

Figure 6. Initial Drive Signals 
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Frequency -Hz 
a) Output PSD Channel 1 
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Frequency -Hz 
b) Output PSD Channel 2 


Figure 10. External Store Response 





